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Abstract 

In a recent publication, we derived the mesoscale continuum theory of plasticity 
for multiple-slip systems of parallel edge dislocations, motivated by the statistical- 
based nonlocal continuum crystal plasticity theory for single-glide due to Yefimov 
et al. (2004b). In this dislocation field theory (DiFT) the transport equations for 
both the total dislocation densities and geometrically necessary dislocation densities 
on each slip system were obtained from the Peach-Koehler interactions through both 
single and pair dislocation correlations. The effect of pair correlation interactions 
manifested itself in the form of a back stress in addition to the external shear and 
the self-consistent internal stress. We here present the study of size effects in single 
crystalline thin films with symmetric double slip using the novel continuum theory. 
Two boundary value problems are analyzed: (1) stress relaxation in thin films on 
substrates subject to thermal loading, and (2) simple shear in constrained films. In 
these problems, earlier discrete dislocation simulations had shown that size effects 
are born out of layers of dislocations developing near constrained interfaces. These 
boundary layers depend on slip orientations and applied loading but are insensitive 
to the film thickness. We investigate stress response to changes in controlled param- 
eters in both problems. Comparisons with previous discrete dislocation simulations 
are discussed. 
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1 Introduction 



Contrary to the prediction of classical crystal plasticity theory, experimental 
observations at length scales ranging from hundreds of nanometers to tens of 
microns show size effects of the type "smaller is harder" (Ebeling and Ashby, 
1966; Brown and Ham, 1971; Fleck et al., 1994; Ma and Clarke, 1995; Stolken 
and Evan, 1998; Arzt, 1998). This failure of conventional continuum theory is 
caused by the lack of a characteristic length scale. Several more sophisticated 
theories (Aifantis, 1984; Walgraef and Aifantis, 1985; Fleck and Hutchinson, 
1993; Fleck et al., 1994; Ortiz and Repctto, 1999; Ortiz et al., 2000; Acharya 
and Bassani, 2000; Acharya and Beaudoin, 2000; Bassini et al., 2001; Gurtin, 
2000, 2002, 2003) have been developed which attempt to incorporate a length 
scale through the concept of geometrically necessary dislocations (GNDs) as 
introduced by Nye (1953). In these theories, however, the length scale enters 
in an ad-hoc fashion, and often has to be supplied a priori by comparison with 
discrete dislocation simulations or experimental results. 

Alternatively, Yefimov et al. (2004b, a) have applied a nonlocal continuum 
plasticity theory based on work by Groma (1997) and Zaiser et al. (2001) 
to successfully solve a set of boundary- value problems for systems with one 
active slip system. ^ They described the evolution of total dislocation densities 
and GND densities using a set of coupled transport equations. In addition to 
external shear and Peach-Koehler interactions among dislocations, the effect 
of pair-dislocation correlation, in the form of a back stress, was considered; 
the latter gave rise to a natural length scale 1 / ^p, determined by the average 
dislocation spacing p. Thriving on the success of their theory, Yefimov and Van 
der Giessen (2005a,b) attempted to extend their single-slip theory to describe 
multiple-slip systems on phenomenological grounds. Albeit favorable results 
were obtained in the problem of shearing of thin films, the theory could not 
capture the size and orientation dependent hardness observed in thin films. 

To address this problem, we have reformulated the multiple-slip theory aiming 
to extract the correct angular dependence of the back stress between different 
pairs of slip orientations (Limkumnerd and Van dcr Giessen, 2007). By solving 
Bogolyubov-Born-Green-Yvon-Kirkwood (BBGYK) integral equations that 
relate different orders of dislocation correlation functions, the functional forms 
of pair-dislocation densities were derived. The results provided slip-orientation 
dependence of pair densities from which the exact expression of the back 
stress was obtained. In their recent publication, Groma et al. (2006) arrived 
at the same expression for a pair correlation function in the case of single-slip 



^ Somewhat similar approaches have been taken by Arsenlis and Parks (2002); 
Arsenlis et al. (2004); El-Azab (2000); Limkumnerd and Sethna (2006); Acharya 
and Roy (2006); Roy and Acharya (2006). 
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systems. 



We begin in Sec. 2 by giving a summary of our continuum theory with a 
short account to the work of Yefimov and Van der Giessen (2005a). In Sec. 3, 
we apply the theory to the problem of stress relaxation in single crystalline 
thin films on substrates subjected to thermal loading. It was this problem in 
which the results between the former multiple-slip theory (Yefimov and Van 
der Giessen, 2005b) and discrete dislocation simulations (Nicola et al., 2003, 
2005b) deviated most. In a quasi-static limit, where dislocations rearrange 
themselves much faster than the stress increase in the film, an analytical solu- 
tion is derived. The hardening effect due to the film thickness and comparisons 
with the discrete dislocation results can be directly investigated for two slip 
orientations. Finally in Sec. 4, we revisit the problem of the simple shear re- 
sponse of thin films, which was used by Yefimov and Van der Giessen (2005a) 
for selecting their slip-interaction law. Layers of dislocations form on the top 
and bottom boundaries which give rise to size effects. Analytical solutions 
of our theory are checked against the discrete dislocation simulations by Shu 
et al. (2001). 



2 Summciry of DiFT-beised plasticity 

Over a decade ago, Groma (1997) has derived a set of transport equations 
governing the motion of many-dislocation densities by carrying out a statistical 
averaging procedure on ensembles of edge dislocations on parallel glide planes. 
Zaiser et al. (2001) later on specialized these equations to describe evolution of 
single-dislocation densities in terms of pair-dislocation densities. Recently the 
authors have extended the above formalism to include systems with more than 
one active slips (Limkumnerd and Van der Giessen, 2007). By constructing 
the integral equations that relate different orders of dislocation correlation 
functions, we explicitly calculate pair correlation functions, and hence pair- 
dislocation densities. In this section we shall briefly summarize this continuum 
theory, leading the derivation to Limkumnerd and Van der Giessen (2007). 

Consider a single crystal with N active slip systems where each system i is 
defined by slip direction Sj and slip plane normal mj. We assume that the 
motion of dislocations is overdamped; positive dislocations on slip system i 
fiow with velocity Vj = (hi/B) rf*^ in the direction of their Burgers vector bj = 
6sj, with magnitude proportional to effective resolved shear stress rf^ with 
drag coefficient B, while negative dislocations flow in the opposite direction. 
The evolution equations for uncorrelated, single-dislocation densities pf and 
P7 can then be re-written in terms of a set of coupled transport equations for 
total dislocation density pi = pf + p^ and the GND density Ki = pf — p^ as 



3 



follows: 



dtPi + V • [niVi] = 0, dtKi + V • [piVi] = 0, (no sum over i) (1) 

with V the derivative with respect to spatial position r. Nucleation and anni- 
hilation of dislocations can be taken into account by modifying the right-hand 
side of the evolution law for pi (cf. Yefimov et al. (2004b)). The dislocation 
density description can be incorporated into the framework of crystal plasticity 
through Orowan's relation 

7. = ^rfA (2) 
and the definition of plastic strain rate 

N 

eP = ^ %Pi , Pi = |(si (8) nij -Firij (g) Sj) . 

i=l 

Substitution into the second dynamical equation in (1) and time integration 
yields Kroner's relation 

ACi = -(l/6)(s,- V)7,, (3) 

which connects GND density Ki to plastic slip ji. 

The effective resolved shear stress 

rf ^n- r^ (4) 

consists of Tj — ^the external shear stress plus the self-consistent, long-range, 
single-dislocation interaction — and the back stress given by 

) - Vcosf^- )^^^-^^^ f5) 

arising from the short-range, dislocation-dislocation interactions. Here fi and 
1/ are the shear modulus and the Poisson ratio respectively. The strength of 
intra-slip back stress is controlled by the dimensionless constant D. The back 
stress contribution from slip system j to slip system i is reduced relative to 
the self back stress by a factor cos(%), where % is the angle between planes 
of slip system i relative to j. 

The form of the back stress as shown in Eq. (5) reduces to that of the previous 
single-shp theory (Groma et al., 2003; Yefimov et al, 2004b,a) for = 1. The 
cos(^jj) slip-interaction coupling considered here also appears in the strain- 
gradient theory for continuum crystal plasticity by Gurtin (Gurtin, 2000, 2002, 
2003). In an early attempt to extend their theory to describe systems with 
multiple slips, Yefimov and Van der Giessen (2005a) had considered three 
different coupling terms: cos^(6'ij), cos(26'ij), and cos(6'jj). They subsequently 
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discarded the first and the third variations upon comparisons with discrete 
dislocation simulations by Shu ct al. (2001). Although the chosen form of 
coupling showed reasonable agreements with the discrete dislocation results 
in the problem of simple shearing of constrained thin film, it failed to capture 
the dependence on film size and slip-orientation observed in the problem of 
stress relaxation in thin films on substrates (Yefimov and Van der Giessen, 
2005b). We shall reexamine these problems with the new continuum theory in 
the following sections and argue that the success of the cos(26'jj)-type coupling 
was just fortuitous. 



3 Application to single crystal thin films on a substrate 

In this section wc consider the problem of stress relaxation in a single crys- 
talline thin film, oriented for symmetric double slip, on a substrate subjected 
to thermal loading. The geometry of the problem is shown in Fig. 1. Initially 
both the thin film, with thermal expansion coefficient ct/, and the substrate, 
with coefficient c^s, are at a (high) temperature Tq. Since a/ > cts, a tensile 
stress up in the film as temperature decreases (with a rate T). At sufficiently 
high stress, pairs of dislocations nucleate on the two active slips according to 
Frank-Read mechanism. When the material is assumed to be initially homo- 
geneous, the problem is effectively one-dimensional; only variations along the 
direction perpendicular to the film matter and the only non- vanishing stress 
component is axx = cr. Also, by symmetry, 71 = —72 and ti = — T2. Hence, 
on average, the density of positive dislocations on the first slip is the same 
as that of negative dislocations on the second, while the negative of the first 
slip and the positive of the second are driven out of the system through the 
top traction-free surface. We shall henceforth drop the subscripts and only 
consider slip system 1. 

This problem can be treated rather simply in a quasi-static limit where disloca- 
tions rearrange themselves much faster than the stress change. In this limit, the 
exact expressions for nucleation and/or annihilation terms are unimportant 
and the nature of the evolution equation (1) is only to transport dislocations 
inside the thin film according to its overall effective stress. At any particular 
time, the distribution of these densities can be calculated from the competition 
between the back stress and the stress due to the thermal mismatch. Given 
the form of the back stress (5), we can derive the time-dependence of this 
expression from the compatibility requirement in terms of slip 7^ on system 1 
and 2. Using Kroner's relation (3), the time evolution of the overall resolved 
shear stress as a function of slip orientation can then be found. The effects of 
film thickness and slip orientation on the stress response can be investigated 
from these expressions. 
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Fig. 1. A thin film of thickness H and thermal expansion coefficient a/ is situated 
on top of an infinite substrate with coefficient ag- The film has two symmetrical 
slip planes defined by angle 0. The y axis is taken to be perpendicular to the 
film-substrate interface. 

In the absence of plasticity, the stress inside the film would build up according 
to 

aN(r)=2//([^)a(To-T), (6) 

where a = af — as \s the effective expansion coefficient of the film relative to 
the substrate. Once the yield point is reached, uy = (7n(7y), plastic straining, 

^xx = -7 sin (20), 

is governed by the resolved shear stress 

sin(2(/)) 



Compatibility of the thermally-induced strain and the elastoplastic strains 
requires that after the yield point, is reached, 

(1 + i/)aAr = ^^p^ (cr - (7y) - 7 sin(20) (7) 

Zi Lb 

where AT = Ty — T is the temperature drop since yield, and 7 is the plastic 
slip (taken to be of shp system 1). The effective shear stress r^^ comprises the 
resolved shear stress 



sin(20) 

T = 



2/i . ,^ , ■ 
o-N + Y3^7sm(20) 



(8) 
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and the back stress which, according to (5), is given by 



27r(l - v)p 
lihDs\n\(t)) 
7r(l — u)k 



sin(0) [1 — cos(20)] dyK 
dyK , 



(9) 



since p — k in this system. Combining eqs. (7)-(9) with k — — sm{(j))/bdyj 
from (3), we can write the effective shear stress as: 



sin(20) 2/i 
2 1-u 



2ii 



V . hD sin3(0) 9^7^ 

crN + 7sm(20) + 



TT sin(20) c?y7 



(10) 



Under the quasi-static assumption mentioned above, the equation of motion 
(1) is solved by force balancing — in other words — by setting r^^ = 0. Eq. (10) 
then gives the nonlinear differential equation 



C/(0) 



7sm(20) + ^^<^N 







with the length scale ( = 2hD/Ti being considered the new fitting parameter 
(instead of D). The solution during yield, subject to the no-slip condition 
7 = at the film-substrate interface = 0, is unique and given by 



sin(20) 7 = 
A = 



1 - z/ 

2/r 

sy 1 



CTn — cry 



aN cosh(A?/) + ay sinh(A|/) 
cry cosh(A|/) -|- (Tn sinh(A|/) 



mc 



(12a) 
(12b) 



Here, 5y = (1 -t- u)a{TQ — Ty) is the film's strain at yield, and /(0) = 
sin^(0)/sin(20) contains the angular dependence on slip orientation. The 
stress profile after yield can be derived using Eq. (8): 



(7(y,r) = (7y 



(Tn cosh (Ay) -|- (jy sinh(A|/) 
(Ty cosh (Ay) -|- (Tn sinh(Ay) 



(13) 



The average stress over the thickness of the film, {uxx) = i^/H) ady, fol- 
lows directly from Eq. (13) as 



try 
XH 



log 



cosh{XH)) + sinh(Aii') 



(Ty 



(14) 



To compare results between the non-local theory and discrete dislocation sim- 
ulations by Nicola et al. (2003) we use parameters from their simulation. The 
film is taken to be isotropic with Poisson's ratio z/ = 0.33, Young modu- 
lus £■ = 70 GPa (from which the value of fj, — E / {2(1 + u)) is computed). 
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and thermal expansion coefficient «/ = 23.2 x 10~^K~^. These values are 
representative of aluminum. The silicon substrate has expansion coefficient 
as = 4.2 X 10~^K~^. The system is cooled from an initial temperature of 
To = 600 K down to T = 400 K, at a rate of T = 4 x 10^ K/s. For the source 
density and source strength (distribution) chosen by Nicola et al. (2003), yield 
starts when the temperature reaches Ty ~ 582 K, i.e. at ay — 35.8 MPa. 

y (/^m) 
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a (MPa) 



Fig. 2. Stress distribution across the film as predicted by Eq. (13) for = 60°. 
Vertical lines indicate the average stress in the films at different film thicknesses H. 

Fitting to the average film stress at T = 400 K predicted by the discrete 
dislocation simulations for orientation = 60° yields a value of C — 28.5 nm. 
Fig. 2 shows the corresponding stress distribution across the film thickness 
according to Eq. (13). At the film-substrate interface, the stress a reaches its 
elastic value of (Tn — 397 MPa, and decays roughly exponentially to the yield 
stress ay ~ 35.8 MPa at the free surface. This profile is independent of the 
film thickness H, as is the discrete dislocation result for the thickest two films. 
The average stress (a^x) for each thickness is indicated by a vertical line. The 
result for = 30° exhibits a similar functional dependence but with a steeper 
decay, and is omitted for brevity. 

Fig. 3 (a) and (b) show the average stress ( function of temperature T 

for different film thicknesses for cf) — 30° and 60° respectively. When the tem- 
perature axis is read right-to-left as a measure of strain, these stress-strain 
curves are steeper (film is harder) as the thickness decreases. The harden- 
ing rate also increases with increasing 0, even though the Schmid factors for 
both orientations are identical. Finally, the prediction of the average tensile 
stress versus film thickness according to Eq. (14) is shown in Fig. 4 against 
the discrete dislocation results (in symbols) for both slip orientations with 
satisfactory agreement. 



The thickness dependence of stress predicted by Eq. (14) is clearly a more 
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H = 0.25 //m 




400 450 500 550 600 (lb) 

Fig. 3. Temperature dependence of average tensile stress {gxx) from Eq. (14) for 
slip orientations (a) = 30° and (b) (f) = 60°. The solid dots represent the discrete 
simulation results (Nicola et al., 2003) 

complicated one than a simple scaling of the type {axx) oc with p varying 
usually between 1/2 and 1. In order to see how large this deviation is, Fig. 5 
shows the data of Fig. 4 on double- log scales. For = 60°, the theoretical 
{(Txx){H) is rather close to a power law over the entire regime considered here, 
but curves upwards for very small H when = 30°. Enhanced hardening in 
very thin films is observed in discrete dislocation results (Nicola et al., 2003, 
2005b) and has been attributed there to dislocation sources being shut down 
by relatively long pile-ups; this effect is absent in the quasi-static solution 
developed here since nucleation is not taken into account. 

A similar theoretical study has been carried out by Nicola et al. (2005a) using 
Gurtin's strain-gradient theory. Compared to the discrete dislocation results. 
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Fig. 4. Average stress at final temperature as a function of film thickness H for 
(p = 30° and 60°. The symbols represent results from the discrete dislocation simu- 
lations. 
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Fig. 5. Log-log plot of the same data as in Fig. 4, and comparison with simple power 
law scaling laws with exponents —0.5 and —1. 

size-dependent hardening was captured but not the orientation dependence 
since Gurtin's original theory predicts the same response for = 30° as for 
= 60°. Subsequently, they proposed a modified "defect energy" based on the 
consideration of dislocation pile-ups which did predict the correct 0-trend. The 
latter implies a material length scale that scales with cos0, while our theory 
predicts scaling with l//(0) oc sin^0/ cos0; the ratio of these for = 30° 
and 60° is identical. It is also interesting to note that the theory by Nicola 
et al. (2005a) reveals a constant hardening rate for a given thickness and shp 
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orientation, whereas we find a weak logarithmic dependence on temperature. 
Both outcomes are within the error bar of the discrete dislocation results. 



4 Simple shecir of constrained film 



We consider the same film as in the previous section, but now subjected to 
a shear F in the x direction, see Fig. 6. While the normal strain was uni- 
form in the film under thermal straining, in the present problem the only 
non-vanishing stress component ct^^j/ is uniform across the width. The sec- 
ond difference is that now both surfaces are impenetrable for dislocations; 
i.e. 7i = 72 = at y = ±H/2 (note that the origin has been placed at the 
center of the film for calculational convenience). By symmetry, ti = T2 and 
7i = 72 which implies that Ki — K2 and pi — P2- We shall therefore omit the 
subscripts. 




Fig. 6. The thin film of thickness H with two impenetrable top and bottom layers 
is under prescribed shearing T(t). The film has two symmetrical slip planes defined 
by angle (f). The origin of the coordinate system is located at the center of the film. 



We can again solve this problem quasi-statically in the manner of Sec. 3. The 
resolved shear stress ti is given by 

T = COS(20) (7^2^, (15) 

while the back stress is, according to Eq. (5), 

= GD/(0) ^ , (16) 

where G = /i6/(27r(l — uj) contains all the material parameters, and /(0) = 
sin(0)(l — cos(20)) — 2sin^(0) captures the slip orientation information. Force 
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balancing, t — — 0, implies that 



Cxy — GD 



f{<t>) dyK 



cos(20) p 



(17) 



Since a^y is uniform across the film thickness by virtue of equilibrium, we 
arrive at the differential equation 



A{(t))dyn{y) = (T:,yp{y) , 



(18) 



above yield. Here A[(j)) = 2GD sm^{4>)/ cos(20) contains the slip orientation 
dependence. 

Under shear, dislocations of one sign (negative when > 7r/4) move towards 
the top y = H/2 where they are blocked, while the opposite-signed dislocations 
will pile-up against the bottom surface; this implies that K.{y) = —sign{y)p{y). 
The solution of Eq. (18) is thus very simple: 



K'iv) = -sign(y)Koexp a^y y/A{(f)) 



(19) 



The constant of integration Kq in general could be a function of the applied 
shear V. Using the relationship (3) between GND density and shp, Eq. (19) 
together with the no-slip boundary conditions give 



7(y) = -7o(r) {1 - exp [X{\y\ - H/2)]} , 



(20) 



where all the integration constants have been absorbed into 7o(r), and 1/A = 
\A{(I))\ / axy(r) gives the approximate characteristic width of the boundary 
layers as a function of the applied shear F. 

Averaging of the decomposition Sij — efj + efj across the width of the sample, 
along with Hooke's law a^y — 2iJ,e^y gives 



axy^2ii (r/2-cos(20) (7)). 



(21) 



Here, we have made use of the fact that {Sxy) =T/2 and employed Eq. (2) to 
find — cos(20)7. The average shp can be calculated directly from Eq. (20), 



(7) = -7o(r) 



-XH/2 



XH/2 



(22) 



The functional form of 70 (E) can be obtained in the limit of large film thick- 
ness, i7 — > 00, where the system is insensitive to the boundary layers which 
results in perfect plasticity. In this case Eq. (21) implies that 



ty = (r + 2 cos(2(/))7o) . 



(23) 
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Eqs. (21)-(23) together provide an implicit expression of Uxy as a function of 
the applied shear F: 




Fig. 7. Comparison of the discrete dislocation (Shu et al., 2001) (dashed lines) 
and nonlocal plasticity (solid curves) shear strain profiles at different values of the 
applied shear V for film thickness H = 1 /xm 

The continuum theory is tested against the discrete dislocation simulations by 
Shu et al. (2001) on a crystal with two shp systems oriented at = 60°. The 
elastic properties are the same as in Sec. 2, i.e. /i = 26.3 GPa and u = 0.33, and 
stress is measured in units of the mean nucleation strength ctq = 1.9 x 10~^/x 
in the discrete simulations. We first note that the width of the boundary 
layers 1/ A cannot be used as a fitting parameter since its value changes with 
increasing stress. We therefore define the length parameter / = o"xy/(o"oA) = 
|A|/(To, which is independent of (T^y, as a new fitting parameter. Given stress 
axy at a selected shear T, the value of / can be determined from fitting Eq. (20) 
to the strain distribution across the film thickness, as shown in Fig. 7. The 
fitting procedure is somewhat intricate due to the non-algebraic nature of 
Eq. (24) which needs to be computed for 'j{y) in Eq. (20) at a given F. We 
therefore take the stress value from the simulation stress-shear curve (Fig. 8) 
as an additional input for the fitting of /, yielding / ~ 46 nm for the case 
of if = 1 iim on the basis of the stress at F = 0.0218. Fig. 7 shows shear 
strain distributions across the film thickness at three other values of F where 
no additional fitting has been performed. 

For further comparison. Fig. 9 (b) shows the theoretical distribution of dislo- 
cation density (recall that p — k ior all slip systems) in comparison with the 
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0.00 0.01 0.02 0.03 0.04 

Fig. 8. Shear stress response to applied shear F for various film thicknesses. The 
discrete dislocation data (dashed lines) is taken from Shu et al. (2001). 

discrete dislocation distribution in a periodic cell with a width of 1 yum. The 
theory correctly predicts the development of intense dislocation boundary lay- 
ers. The core of the crystal is left almost dislocation free as dislocations pass 
each other almost unhindered on their way towards the top and bottom faces. 



From the above-mentioned best-fit / for the film thickness of if = 1 /im, we can 
study the shear response for different film thicknesses. Data of the discrete dis- 
location simulations suggest thickness-dependent initial yield strengths. The 




Fig. 9. Dislocation distribution in a = 1 /im thick (j) = 60° film at an overall shear 
of F = 0.0218 according to (a) discrete simulations by Shu et al. (2001) and (b) the 
theoretical solution (19) with / ~ 46 nm. 
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Fig. 10. Dislocation distribution in a i7 = 1 /xm thick (f) = 30° film at an overall 
shear of F = 0.0218 according to (a) discrete simulations by Shu et al. (2001) and 
(b) the theoretical solution (19) with / ~ 8.7 nm. 

responses are shown in Fig. 8 in comparison with results from the discrete sim- 
ulations. We supply for each film thickness the best-fit yield point as an extra 

degree of freedom. Similar to the previous test problem (Sec. 3), the stress- 
strain curves show size-dcpcndcnt hardening. The hardening rate decreases 
with increasing applied external shear, and approaches a constant value at 
large shear. Shu et al. (2001) also analyzed this problem with their version 
of strain-gradient theory and found weak size effects. Their stress response, 
however, is linear due to the fact that the width of dislocation boundary lay- 
ers is constant in their theory. The same linear stress-strain relation was also 
predicted by Gurtin's strain-gradient theory (Bittencourt et al., 2003). 

It should be mentioned that the exact form of the slip-interaction coupling 

(in Eq. (5)) turns out to be unimportant in this problem. The slip orientation 
dependence is buried inside the definition of |A(0)| which has been absorbed 
into the fitting parameter /. On this ground, it does not matter whether this 
coupling be cos(^y) or cos(26'jj) as proposed by Yefimov and Van der Giessen 
(2005a). 

Our theory predicts drastic changes in behavior when crosses 45°. Due to a 
sign change in the resolved shear stress, the charges of dislocations at the two 
interfaces reverse from the present situation when < 45° (resulting in the 
sign alternations of Kq and 70 in Eqs. (19), (20), (22), and (23)). As a result, 
the applied shear acts in favor of the new dislocation arrangement — in other 
words — our theory predicts that the back stress further enhances plasticity 
instead of impeding the fiow of new dislocations into the boundaries. Hence, 
thinner boundary layers are expected which suggests smaller size effects. More 
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Fig. 11. Types of dislocations at the pile-ups for (a) (f) = 30° and (b) (f) = 60°, and 
their locking situations as they glide (a) upwards or downwards (as viewed upwide 
down) for 30° case or (d) in the opposite directions for 60° case. 

quantitatively, for the orientation angle of, say, = 30°, the layers should be 
thinner by a factor of ho^/ko^ = |A(30°)/A(60°)| = (sin(30°)/ sin(60°))3 ~ 
0.19. The dislocation distribution thus predicted is shown in Fig. 10(b). 

Discrete dislocation dynamics simulations, however, reveal essentially no bound- 
ary layers at all — or, equivalently, boundary layers that span the entire width 
(Fig. 10(a)). Upon closer examination, we find 'locks' of dislocations^ of like 
charges on different slip systems which prevent their motion pile-ups to the 
boundaries. A pair of dislocations with the relative angle of their Burgers 
vectors between 90° and 270° feel their mutual attraction when they ghde 
past each other. Although rather weak, this interaction is apparently strong 
enough in this case for locking to occur. Figs. 11(a) and (b) show the types 
of dislocations which accumulate at the boundaries for = 30° and 60°, re- 
spectively. Figs. 11(c) and (d) demonstrate one of the two situations when 
locking happens in each case (the others are 180° rotations of these). In the 
region sufficiently far away from the boundaries, event 11(c) is roughly as 
likely to occur as event 11(d), since the situations differ just by a 90° rotation 

^ The use of the phrase 'locks' for parallel edge dislocations is somewhat inappro- 
priate, because such dislocations only interact through their long-range field but do 
not alter the topology as happens in, e.g., Lomer locks. We nevertheless use the 
term here because it expresses the small-scale interactions that obstruct dislocation 
motion. 
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followed by a flip about the y-axis. The relative likelihood, however, increases 
immensely close to the boundaries because only in the 30° case do dislocations 
moving to the same boundary permit locking, Fig. 11(c). This mutual locking 
of slip systems prevents dislocations to reach the boundaries and form local- 
ized boundary layers. The locking mechanism is purely a discrete phenomena 
and cannot be captured by the current continuum theory without further re- 
finement. Due to its relatively small probability, locking seldom occurs in the 
60° case. 



5 Discussion and Conclusion 

We applied the recently formulated multislip continuum plasticity theory to 
analyze two boundary value problems relating to thin films. In Sec. 3, we 
studied stress relaxation mechanism in thin films on substrates with thermal 
loading. Wc obtained an explicit analytical expression of the stress distribution 
as a function of slip orientation with one fitting parameter. The predictions 
were in good agreement with the discrete dislocation results of Nicola et al. 
(2003, 2005b). Our theory was able to show size-dependent hardening and 
the hardening due to slip orientations — both of which the previous continuum 
theory failed to explain. Subsequently, we analyzed simple shear in constrained 
films. Similarly to the first problem, we observed dislocation pile-ups at the 
top and bottom constrained boundaries. The thickness of dislocation layers 
depends weakly on the incremental shear. Regardless of the difference between 
the forms of slip-interaction coupling between our theory and that in Yefimov 
and Van der Giessen (2005a), our theory also gave satisfactory agreements 
with results from discrete dislocation dynamics simulations (Shu et al., 2001). 
We pointed out that this term can be absorbed into fitting parameter; the 
correct functional form of the coupling, therefore, cannot be decided only on 
the basis of this problem. 
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